home *** CD-ROM | disk | FTP | other *** search
/ NeXT Education Software Sampler 1992 Fall / NeXT Education Software Sampler 1992 Fall.iso / Mathematics / Notebooks / SigProc2.0 / Packages / SignalProcessing / ObjectOriented / RewriteRules.m < prev    next >
Encoding:
Text File  |  1992-08-18  |  20.2 KB  |  637 lines
  1. (*  :Title:    Rewrite Rules for Signal Processing Expressions  *)
  2.  
  3. (*  :Authors:    Brian Evans, James McClellan  *)
  4.  
  5. (*  :Summary:    To provide system rewrite rules a la Myers and Covell.  *)
  6.  
  7. (*  :Context:    SignalProcessing`ObjectOriented`RewriteRules`  *)
  8.  
  9. (*  :PackageVersion:  2.7    *)
  10.  
  11. (*
  12.     :Copyright:    Copyright 1989-1992 by Brian L. Evans
  13.         Georgia Tech Research Corporation
  14.  
  15.     Permission to use, copy, modify, and distribute this software
  16.     and its documentation for any purpose and without fee is
  17.     hereby granted, provided that the above copyright notice
  18.     appear in all copies and that both that copyright notice and
  19.     this permission notice appear in supporting documentation,
  20.     and that the name of the Georgia Tech Research Corporation,
  21.     Georgia Tech, or Georgia Institute of Technology not be used
  22.     in advertising or publicity pertaining to distribution of the
  23.     software without specific, written prior permission.  Georgia
  24.     Tech makes no representations about the suitability of this
  25.     software for any purpose.  It is provided "as is" without
  26.     express or implied warranty.
  27.  *)
  28.  
  29. (*  :History:    *)
  30.  
  31. (*  :Keywords:    *)
  32.  
  33. (*
  34.     :Source:    R. H. Bamberger.  {The Directional Filter Bank:
  35.           A Multirate Filter Bank for the Directional Decomposition
  36.           of Images},  Georgia Tech Ph. D. Thesis.  November 1990.
  37.  
  38.         M. M. Covell.  {An Algorithm Design Environment for
  39.           Signal Processing}.  M.I.T. Ph. D. Thesis.  December,
  40.           1989.
  41.  
  42.         R. Crochiere and L. Rabiner, {Multirate Digital Signal
  43.           Processing}.  Prentice Hall: Englewood Cliffs, NJ.  1983.
  44.  
  45.         B. L. Evans, R. H. Bamberger, and J. H. McClellan,
  46.           "Rules for Multidimensional Multirate Structures,"
  47.           to be published in the {IEEE Journal on Signal Processing}.
  48.  
  49.         J. Kovacevic and M. Veterli, "Perfect Reconstruction Filter
  50.           Banks with Rational Sampling Rates in One- and
  51.           Two-dimensions," in Proc. SPIE Conference on Visual
  52.           Communications and Image Processing (Philadelphia, PA),
  53.           pp. 1258-1268, November 1989.
  54.  
  55.         C. S. Myers.  {Signal Representation for Symbolic and
  56.           Numeric Processing}.  M.I.T. Ph. D. Thesis.  August,
  57.           1986.  Appendix D.
  58.  
  59.         P. P. Vaidyanthan and T. Chen, "Some Fundamental Issues
  60.           in Multidimensional Multirate Signal Processing,"
  61.           to be published in {IEEE Journal on Signal Processing}.
  62.  *)
  63.  
  64. (*  :Warning:    *)
  65.  
  66. (*  :Mathematica Version:  1.2 or 2.0  *)
  67.  
  68. (*  :Limitation:  *)
  69.  
  70. (*
  71.     :Discussion:  Systems are defined in the same format as are mathematical
  72.           operators.  This format was chosen so that Mathematica's
  73.           powerful pattern matcher could be used without adjustment.
  74.           A full blown object representation was not needed since
  75.           no cacheing, fetching, dynamic properties, etc. are needed
  76.           as they are in signals.  That is, a z-transform operator
  77.           always exhibits linearity, and its properties do not
  78.           change.  In addition, systems are not instantiated as
  79.           are objects.
  80.  
  81.           Properties are attached to the head of a system name.
  82.           The system rewrite rules need these properties to
  83.           rewrite expressions.
  84.  
  85.     Imag[x_ y_] :> Real[x] Imag[y] + Imag[x] Real[y]
  86.     Real[x_ y_] :> Real[x] Real[y] - Imag[x] Imag[y]
  87.  *)
  88.  
  89. (*  :Functions:     *)
  90.  
  91.  
  92.  
  93. (*  B E G I N     P A C K A G E  *)
  94.  
  95. BeginPackage[ "SignalProcessing`ObjectOriented`RewriteRules`",
  96.           "SignalProcessing`ObjectOriented`System`",
  97.           "SignalProcessing`Support`FilterSupport`",
  98.           "SignalProcessing`Support`LatticeTheory`",
  99.           "SignalProcessing`Support`SigProc`",
  100.           "SignalProcessing`Support`SupCode`" ]
  101.  
  102.  
  103. If [ TrueQ[ $VersionNumber >= 2.0 ],
  104.      Off[ General::spell ];
  105.      Off[ General::spell1 ] ];
  106.  
  107.  
  108. (*  U S A G E     I N F O R M A T I O N  *)
  109.  
  110. SPRecursiveRewrite::usage =
  111.     "SPRecursiveRewrite[expr] only returns a fully rewritten version of \
  112.     the signal processing expression expr. \
  113.     The Rewrite knowledge base is let loose on the expression expr \
  114.     without any guidance. \
  115.     The user can see the intermediate expressions by setting the \
  116.     Dialogue option to True or False."
  117.  
  118. SystemRewriteRules::usage =
  119.     "SystemRewriteRules contains a list of over 50 rewrite rules for \
  120.     operators (systems) . \
  121.     This rule is composed of three other rule bases: \
  122.     (1) rules based on operator properties, \
  123.     (2) rules for single rate digital signal processing, and \
  124.     (3) rules for multirate digital signal processing. \
  125.     The rules in (1) are patterned after those found in Myers' Signal \
  126.     Processing Language and Interactive Computing Environment (SPLICE) \
  127.     and Covell's Algorithm Design Environment (ADE)."
  128.  
  129. (*  E N D     U S A G E     I N F O R M A T I O N  *)
  130.  
  131.  
  132. Begin[ "`Private`" ]
  133.  
  134.  
  135. (* commutativemult --  multiplication always commutes in 1-D *)
  136. commutativemult[x_, y_, n_Symbol] := True
  137. commutativemult[x_, y_, n_List] := ( matmult[x, y, n] == matmult[y, x, n] )
  138.  
  139. (* extractit -- undoes nesting operation *)
  140. extractit[ a_, b__, h_ ] := ( extractit[a, h], extractit[b, h] )
  141. extractit[ h_[a__], h_ ] := ( extractit[a, h] )
  142. extractit[ h_[a__], hother_ ] := h[a] /; ! SameQ[h, hother]
  143.  
  144. (* gcd  *)
  145. gcd[{m_Integer, n_Integer}] := GCD[m, n]
  146. gcd[m_Integer, n_Integer] := GCD[m, n]
  147.  
  148. (* getdiagonal *)
  149. getdiagonal[mat_] := getdiagonal[mat, Min[Dimensions[mat]]]
  150. getdiagonal[mat_, n_] :=
  151.     Block [    {i},
  152.         Table[ mat[[i,i]], {i, 1, n} ] ]
  153.  
  154. (* isDownsamplerQ *)
  155. isDownsamplerQ[ PolyphaseDownsample ] = True
  156. isDownsamplerQ[ Downsample ] = True
  157. isDownsamplerQ[ x_ ] := False
  158.  
  159. (* isLSIfilterQ *)
  160. isLSIfilterQ[ FIR ] = True
  161. isLSIfilterQ[ IIR ] = True
  162. isLSIfilterQ[ x_ ] := False
  163.  
  164. (* isUpsamplerQ *)
  165. isUpsamplerQ[ PolyphaseUpsample ] = True
  166. isUpsamplerQ[ Upsample ] = True
  167. isUpsamplerQ[ x_ ] := False
  168.  
  169. (* iszero  *)
  170. iszero[ l_Symbol ] := SameQ[l, 0]
  171. iszero[ l_List ] := SameQ[l, Apply[Table, {0, {Length[l]}}]]
  172.  
  173. (* lsi *)
  174. lsi[ x_ ] := LINEAR[x] && SHIFTINVARIANT[x]
  175.  
  176. (* matmult *)
  177. matmult[a_, b_, n_List] := a . b
  178. matmult[a_, b_, n_Symbol] := a b
  179.  
  180. (* matinv  *)
  181. matinv[a_, n_List] := Inverse[a]
  182. matinv[a_, n_Symbol] := 1/a
  183.  
  184. (* mod *)
  185. mod[l_Integer, L_Integer] := Mod[l, L]
  186. mod[l_?IntegerVectorQ, L_?ResamplingMatrixQ] := ResamplingMatrixMod[l, L]
  187.  
  188. (* modulatorq *)
  189. modulatorq[f_, n_Symbol] := ! FreeQ[f, n]
  190. modulatorq[f_, {n_Symbol}] := ! FreeQ[f, n]
  191. modulatorq[f_, {n_Symbol, nrest__}] :=
  192.     (! FreeQ[f, n]) || modulatorq[f, {nrest}]
  193.  
  194. (* removecf --  remove common factors in separable or 1-D resampling *)
  195. removecf[m1_, m2_, n_Symbol] :=
  196.     Block [    {gcdvalue},
  197.         gcdvalue = gcd[m1, m2];
  198.         {m1 / gcdvalue, m2 / gcdvalue} ]
  199.  
  200. removecf[diagmat1_, diagmat2_, n_List] :=
  201.     Block [    {diag1, diag2, gcdlist},
  202.         diag1 = getdiagonal[diagmat1];
  203.         diag2 = getdiagonal[diagmat2];
  204.         gcdlist = Map[gcd, Transpose[{diag1, diag2}]];
  205.         { DiagonalMatrix[diag1 / gcdlist],
  206.           DiagonalMatrix[diag2 / gcdlist] } ]
  207.  
  208. (* resampledfilter *)
  209. resampledfilter[head_, factor_Integer, n_Symbol, coeffs_List, rest___] :=
  210.     Block [    {i, numcoeffs},
  211.         numcoeffs = Length[coeffs];
  212.         newcoeffs = Table[ coeffs[[i]], {i, 1, numcoeffs, factor} ];
  213.         head[n, newcoeffs, rest] ]
  214.  
  215. (* resampledcoeffs *)
  216. resampledcoeffs[factor_Integer, coeffs_List] :=
  217.     Block [    {i, matchlist, numcoeffs},
  218.         numcoeffs = Length[coeffs];
  219.         matchlist = Table[0, {numcoeffs}];
  220.         For [ i = 1, i <= numcoeffs, i += factor, matchlist[[i]] = _ ];
  221.         MatchQ[coeffs, matchlist] ]
  222.  
  223. (* sepresq -- detect separable or 1-D resampling *)
  224. sepresq[mat_, n_Symbol] := IntegerQ[mat]
  225. sepresq[mat_, n_List] := DiagonalMatrixQ[mat] && ResamplingMatrix[mat]
  226.  
  227. (* shiftq *)
  228. shiftq[expr_] := MatchQ[expr, Shift[l_, n_][x_]]
  229.  
  230. (* switchtest *)
  231. switchtest[L_, M_, Right] :=
  232.     Block [    {cond, dL, dM, uLinv, uMinv, vLinv, vMinv},
  233.         {cond, {uMinv, dM, vMinv}, {uLinv, dL, vLinv}} =
  234.           SmithFormSameV[SmithNormalForm, M, L];
  235.         cond && RelativelyPrime[dL, dM] ] /;
  236.     ResamplingMatrix[L] && ResamplingMatrix[M]
  237.  
  238. switchtest[L_, M_, Left] :=
  239.     First[ SmithFormSameU[SmithNormalForm, M, L] ] /;
  240.     ResamplingMatrix[L] && ResamplingMatrix[M]
  241.  
  242.  
  243. (*  G L O B A L S  *)
  244.  
  245. RulesBasedOnProperties = {
  246.  
  247.     (*  Shift invariant operator applied to shift  [Myers, 216]    *)
  248.     t_ [ Shift[l_, n_][x_] ] :>
  249.         Shift[l,n] [ t [ x ] ] /; SHIFTINVARIANT[t] ,
  250.     t_ [ s1_, s__ ] :>
  251.         Block [    {signals, signalstomove},
  252.             signals = {s1, s};
  253.             signalstomove = Select[signals, shiftq];
  254.             MoveOperatorsToFront[t, signals,
  255.                          signalstomove, Shift] ] /;
  256.         SHIFTINVARIANT[t] && Count[{s1, s}, Shift[l_,n_][x_]] > 0 ,
  257.  
  258.     (*  Memoryless operator applied to shifts  [Myers, 216]        *)
  259.     (*    t_ [ Shift[l_, n_] [x_] ] :> Shift[l, n] [ t [ x ] ] /;    *)
  260.     (*    MEMORYLESS[t],                        *)
  261.     (*  is handled by the one-input to a shift-invariant system    *)
  262.     (*  rule above because one-input memoryless systems are SI.    *) 
  263.     t_ [ Shift[l_, n_] [x_], y__ ] :> Shift[l, n] [ t [ x, GetArgs[y] ] ] /;
  264.         MEMORYLESS[t] && SameFormQ[Shift[l,n][s_], y] ,
  265.  
  266.     (*  Memoryless operator applied to downsample  [Myers, 216]    *)
  267.     t_ [ Downsample[m_, n_] [ x_ ] ] :> Downsample[m,n] [ t[x] ] /;
  268.         MEMORYLESS[t] ,
  269.     t_ [ Downsample[m_, n_] [ x_ ], y__ ] :>
  270.         Downsample[m,n] [ t [ x, GetArgs[y] ] ] /;
  271.         MEMORYLESS[t] && SameFormQ[Downsample[m, n][s__], y] ,
  272.  
  273.     (*  Memoryless operator applied to upsample  [Myers, 216]    *)
  274.     t_ [ Upsample[l_, n_] [ x_ ] ] :>
  275.         Upsample[l, n] [ t[ x ] ] /;
  276.         MEMORYLESS[t] ,
  277.     t_ [ Upsample[l_, n_] [ x_ ], y__ ] :>
  278.         Upsample[l, n] [ t [ x, GetArgs[y] ] ] /;
  279.         MEMORYLESS[t] && SameFormQ[Upsample[l, n][s__], y] ,
  280.  
  281.     (*  Additive system and sum input  [Myers, 217]        *)
  282.     t_ [ x_ + y_ ] :>  t[x] + t[y] /; ADDITIVE[t] ,
  283.  
  284.     (*  Additive system and zero input  [Myers, 217]    *)
  285.     t_ [ zero_ ] :> 0 /; ADDITIVE[t] && ZeroQ[zero] ,
  286.  
  287.     (*  Homogeneous system and scaled input  [Myers, 217]    *)
  288.     t_ [ g_ x_. ] :> g t[x] /;
  289.         HOMOGENEOUS[t] && (! SameQ[g, x, 1]) &&
  290.         Implies[ ! AtomQ[t],
  291.              MyFreeQ[g, GetOperatorVariables[t]] ] ,
  292.  
  293.     (*  Commutative system reorder  [Myers, 216]    *)
  294.     t_ [ x_, y_] :> t [ y, x ] /; COMMUTATIVE[t] ,
  295.  
  296.     (*  Associative system reorder  [Myers, 217]    *)
  297.     t_ [ x_, t_ [ y_, w_ ] ] :> t [ t [ x, y ], w ] /; ASSOCIATIVE[t] ,
  298.  
  299.     (*  Linear shift-invariant props not in CM's thesis [O & W, 85-86]  *)
  300.     t_ [ u_ [ x_ ] ] :> u [ t[x] ] /; lsi[t] && lsi[u]
  301.       
  302. }
  303.  
  304.  
  305. SingleRateRules = {
  306.  
  307.     (*  Collect cascaded convolution operations under one operator    *)
  308.     Convolve[n_] [ x___, Convolve[n_][y1_, yrest__], z___  ] :>
  309.         Convolve[n] [ extractit[ x, y1, yrest, z, Convolve[n] ] ] ,
  310.  
  311.     (*  Convolution of upsampled signals [Myers, 221]  *)
  312.     Convolve[n_] [ Upsample[l_, n_][x_], y__ ] :>
  313.         Upsample[l, n] [ Convolve[n] [ x, GetArgs[y] ] ] /;
  314.         SameFormQ[ Upsample[l,n][s__], y ] ,
  315.  
  316.     (*  Combining redundant operations [Myers, 217/220] [Bamberge]  *)
  317.     Periodic[m1_, w_] [ Periodic[m2_, w_] [x_] ] :>
  318.         Periodic[m1 m2, w][x] ,
  319.     ScaleAxis[l1_, w_] [ ScaleAxis[l2_, w_] [x_] ] :>
  320.         ScaleAxis[l1 l2, w][x] ,
  321.     Shift[m1_, n_] [ Shift[m2_, n_] [x_] ] :>
  322.         Shift[m1 + m2, n][x] ,
  323.  
  324.     (*  1-D Shift of a periodic signal-- could repeat ad infinitum  *)
  325.     Shift[m_, n_Symbol] [ Periodic[k_, n_Symbol] [x_] ] :>
  326.         Shift[ Mod[m,k], n] [ Periodic [k,n] [x] ] /;
  327.         ( ! SameQ[ m, Mod[m, k] ] ) && ( ! MatchQ[ m, Mod[f_, k] ] ) ,
  328.  
  329.     (*  Conjugate of a product [Myers, 218]  *)
  330.     Conjugate[ x_ y_ ] :> Conjugate[x] Conjugate[y] ,
  331.  
  332.     (*  Real part of a product [Myers, 218] *)
  333.     Re[ x_ y_ ] :> Re[x] Re[y] - Im[x] Im[y] ,
  334.  
  335.     (*  Imaginary part of a product [Myers, 218] *)
  336.     Im[ x_ y_ ] :> Re[x] Im[y] + Im[x] Re[y] ,
  337.  
  338.     (*  Extended math functions [Myers, 220]  *)
  339.     (*  signal + 0 = signal                built-in    *)
  340.     (*  a^n1 a^n2 = a^(n1 + n2)            built-in    *)
  341.  
  342.     (*  Factorizations [Myers, 221] *)
  343.     x_ + g_ x_ :> (1 + g) x ,
  344.     g1_ x_ + g2_ x_ :> (g1 + g2) x    (* also covers g x + g y = g (x + y) *)
  345. }
  346.  
  347.  
  348. MultirateRules = {
  349.  
  350.     (*  Identities for downsamplers [Myers, 217] [Evans, fig. 6] *)
  351.  
  352.     Shift[n0_, n_] [ Downsample[M_, n_] [x_] ] :>
  353.         Downsample[M, n] [ Shift[matmult[M, n0, n], n][x] ] ,
  354.  
  355.     Downsample[M_, n_] [ Shift[n0_, n_] [x_] ] :>
  356.         Shift[matmult[matinv[M, n], n0, n], n] [ Downsample[M, n][x] ] ,
  357.  
  358.     phi_ Downsample[M_, n_][x_] :>
  359.         Downsample[M, n][x Upsample[M, n][phi]] /;
  360.         modulatorq[phi, n] ,
  361.  
  362.     Downsample[M_, n_][ phi_ x_ ] :>
  363.         Downsample[M, n][x] Downsample[M, n][phi] /;
  364.         modulatorq[phi, n] ,
  365.  
  366.  
  367.     (*  Identities for upsamplers [Myers, 219] [Evans, fig. 7] *)
  368.  
  369.     Upsample[L_, n_] [ Shift[n0_, n_][x_] ] :>
  370.         Shift[matmult[L, n0, n], n][ Upsample[L, n][x] ] ,
  371.  
  372.     Shift[n0_, n_][ Upsample[L_, n_] [ x_] ] :>
  373.         Upsample[L, n][ Shift[matmult[matinv[L, n], n0, n], n][x] ] ,
  374.  
  375.     phi_ Upsample[L_, n_][x_] :>
  376.         Upsample[L, n][x Downsample[L, n][phi]] /;
  377.         modulatorq[phi, n] ,
  378.  
  379.     Upsample[L_, n_][ phi_ x_ ] :>
  380.         Upsample[L, n][phi] Upsample[L, n][x] /;
  381.         modulatorq[phi, n] ,
  382.  
  383.  
  384.     (*  Identities for cascades of up and downsamplers [Evans, fig. 9] *)
  385.  
  386.     Downsample[S_, n_] [ Upsample[S_, n_][x_] ] :> x ,
  387.  
  388.     Upsample[S_, n_] [ Downsample[S_, n_][x_] ] :>
  389.         x ImpulseTrain[S, n] ,
  390.  
  391.     Downsample[M1_, n_] [ Downsample[M2_, n_] [x_] ] :>
  392.         Downsample[matmult[M1, M2, n], n][x] ,
  393.  
  394.     Upsample[L1_, n_] [ Upsample[L2_, n_] [x_] ] :>
  395.         Upsample[matmult[L1, L2, n], n][x] ,
  396.  
  397.     Downsample[M_Integer, n_Symbol] [ Upsample[L_Integer, n_Symbol][x_] ] :>
  398.         Upsample[ L / M, n ][x] /;
  399.         IntegerQ[ L / M ],
  400.  
  401.     Downsample[M_List, n_List] [ Upsample[L_List, n_List][x_] ] :>
  402.         Upsample[Inverse[M] . L, n][x] /;
  403.         IntegerVectorQ[ Inverse[M] . L ] ,
  404.  
  405.  
  406.     (*  Commutativity of cascades of up and downsamplers [Evans, fig. 10] *)
  407.     (*  Phase out last two rules when RelativelyPrime is complete  *)
  408.  
  409.     Downsample[M1_, n_] [ Downsample[M2_, n_][x_] ] :>
  410.         Downsample[M2, n] [ Downsample[M1, n][x] ] /;
  411.         commutativemult[M1, M2, n] ,  (* always commutes in 1-D *)
  412.  
  413.     Upsample[L1_, n_] [ Upsample[L2_, n_][x_] ] :>
  414.         Upsample[L2, n] [ Upsample[L1, n][x] ] /;
  415.         commutativemult[L1, L2, n] ,  (* always commutes in 1-D *)
  416.  
  417.     Downsample[M_, n_] [ Upsample[L_, n_][x_] ] :>
  418.         Upsample[L, n] [ Downsample[M, n] [x] ] /;
  419.         commutativemult[M, L, n] && RelativelyPrime[L, M] ,
  420.  
  421.     Upsample[L_, n_] [ Downsample[M_, n_][x_] ] :>
  422.         Downsample[M, n] [ Upsample[L, n] [x] ] /;
  423.         commutativemult[L, M, n] && RelativelyPrime[L, M] ,
  424.  
  425.     Downsample[M_, n_List] [ Upsample[L_, n_][x_] ] :>
  426.         Upsample[L, n] [ Downsample[M, n] [x] ] /;
  427.         CommutableResamplingMatricesQ[M, L] ,
  428.  
  429.     Upsample[L_, n_List] [ Downsample[M_, n_][x_] ] :>
  430.         Downsample[M, n] [ Upsample[L, n] [x] ] /;
  431.         CommutableResamplingMatricesQ[L, M] ,
  432.  
  433.  
  434.     (*  Fundamental Identities Based on Smith Form Dec. [Evans, fig. 11] *)
  435.  
  436.     Upsample[U_, n_List][x_] :> Downsample[Inverse[U], n][x] /;
  437.         RegUnimodularResMatrixQ[U] ,
  438.  
  439.     Downsample[U_, n_List][x_] :> Upsample[Inverse[U], n][x] /;
  440.         RegUnimodularResMatrixQ[U] ,
  441.  
  442.     Downsample[V_, n_List] [ Upsample[U_, n_][x_] ] :>
  443.         Upsample[Inverse[V] . U, n] [x] /;
  444.         RegUnimodularResMatrixQ[U] && RegUnimodularResMatrixQ[V] ,
  445.  
  446.     Upsample[U_, n_List] [ Downsample[V_, n_][x_] ] :>
  447.         Downsample[V . Inverse[U], n] [x] /;
  448.         RegUnimodularResMatrixQ[U] && RegUnimodularResMatrixQ[V] ,
  449.  
  450.  
  451.     (*  Removing redundancy in up/downsamplers cascades [Evans, fig. 12]  *)
  452.  
  453.     Downsample[M_, n_] [ Upsample[L_, n_][x_] ] :>
  454.         Block [    {newL, newM},
  455.             {newL, newM} = removecf[L, M, n];
  456.             Downsample[newM, n][ Upsample[newL, n][x] ] ] /;
  457.         sepresq[M, n] && sepresq[L, n] && ! RelativelyPrime[M, L] ,
  458.  
  459.     Upsample[L_, n_List] [ Downsample[M_, n_List][x_] ] :>
  460.         Block [    {cond, dL, dM, uLinv, uMinv, vLinv, vMinv},
  461.             {cond, {uMinv, dM, vMinv}, {uLinv, dL, vLinv}} =
  462.               SmithFormSameV[SmithNormalForm, M, L];
  463.             newM = Inverse[uMinv] . dM;
  464.             newL = Inverse[uLinv] . dL;
  465.             Upsample[newL, n][ Downsample[newM, n][x] ] ] /;
  466.         First[ SmithFormSameV[SmithNormalForm, M, L] ] ,
  467.  
  468.     Downsample[M_, n_List] [ Upsample[L_, n_List][x_] ] :>
  469.         Block [    {cond, dL, dM, uLinv, uMinv, vLinv, vMinv},
  470.             {cond, {uLinv, dL, vLinv}, {uMinv, dM, vMinv}} =
  471.               SmithFormSameV[SmithNormalForm, L, M];
  472.             newL = Inverse[uLinv] . dL;
  473.             newM = Inverse[uMinv] . dM;
  474.             Downsample[newM, n][ Upsample[newL, n][x] ] ] /;
  475.         First[ SmithFormSameU[SmithNormalForm, L, M] ] ,
  476.  
  477.  
  478.     (*  Switch operations in a non-commutable cascade [Evans, fig. 13] *)
  479.  
  480.     Upsample[L_, n_List][ Downsample[M_, n_List][x_] ] :>
  481.         Block [    {cond, dL, dM, uLinv, uMinv, vLinv, vMinv},
  482.             {cond, {uMinv, dM, vMinv}, {uLinv, dL, vLinv}} =
  483.               SmithFormSameV[SmithNormalForm, M, L];
  484.             newM = dM . uLinv;
  485.             newL = dL . uMinv;
  486.             Downsample[newM, n][ Upsample[newL, n][x] ] ] /;
  487.         switchtest[L, M, Right] ,
  488.  
  489.     Downsample[M_, n_List][ Upsample[L_, n_List][x_] ] :>
  490.         Block [    {cond, dim, dL, dM, i, uLinv, uMinv, vLinv, vMinv},
  491.             {cond, {uMinv, dM, vMinv}, {uLinv, dL, vLinv}} =
  492.               SmithFormSameU[SmithNormalForm, M, L];
  493.             dim = Length[n];
  494.             For [ i = 1, i <= dim, i++,
  495.                   dcommon = GCD[ dM[[i,i]], dL[[i,i]] ];
  496.                   dM[[i,i]] /= dcommon;
  497.                   dL[[i,i]] /= dcommon ];
  498.             newM = vLinv . dM;
  499.             newL = vMinv . dL;
  500.             Upsample[newL, n][ Downsample[newM, n][x] ] ] /;
  501.         switchtest[M, L, Left] ,
  502.  
  503.     (*  Downsampled of shifted upsample [Myers, 219]  *)
  504.     Downsample[L_, n_] [ Shift[l_, n_] [ Upsample[L_, n_][x_] ] ] :>
  505.         Block [ {n0},
  506.             n0 = matmul[matinv[L, n], l, n];
  507.             Shift[n0, n][x] ] /;
  508.         iszero[ mod[l, L] ] ,
  509.  
  510.     Downsample[L_, n_] [ Shift[l_, n_] [ Upsample[L_, n_][x_] ] ] :>
  511.         0 /;
  512.         ! iszero[ mod[l, L] ] ,
  513.  
  514.     (* Downsample of shifted upsample with co-prime M and L [Myers, 219] *)
  515.     Downsample[M_, n_] [ Shift[l_, n_] [ Upsample[L_, n_][x_] ] ] :>
  516.         Block [    {l1, l2},
  517.             {l1, l2} = EuclidFactors[l, M, L];
  518.             Shift[l1, n] [
  519.               Upsample[L, n] [
  520.                 Downsample[M, n] [
  521.                   Shift[l2, n] [ x ] ] ] ] ] /;
  522.         RelativelyPrime[M, L, Right] ,
  523.  
  524.     (*  Interleave as sum of upsampled signals  [Myers, 219]  *)
  525.     Summation[i_, 1, L_, 1] [
  526.       Shift[i_, n_] [
  527.         Upsample[L_, n_] [
  528.           Element[x_List, i_] ] ] ] :>
  529.         Interleave[n][ GetArgs[ToList[x]] ] /;
  530.         SameQ[ Length[x], L ] ,
  531.  
  532.  
  533.     (*  1-D Downsample of interleave (L rel. prime to m)  [Myers, 220]  *)
  534.     Downsample[m_, n_Symbol] [Interleave[n_Symbol][x1_, x__]] :>
  535.         Block [    {ii, len, xlist},
  536.             xlist = List[x1, x];
  537.             len = Length[xlist];
  538.             oneterm = Downsample[m, n] [
  539.                     Shift[- ii m / len, n] [
  540.                       Element[ xlist, Mod[ii m, len] + 1 ] ] ];
  541.             Apply[Interleave[n], Table[oneterm, {ii,0,len-1}]] ] /;
  542.         gcd[Length[{x1, x}], m] == 1 ,
  543.  
  544.     (*  Switching filter and decimator/interpolator [C & R]  *)
  545.     (*  -- pattern matching assumes that the Downsample and  *)
  546.     (*     and PolyphaseDownsample takes the same first and  *)
  547.     (*     last parameter, but the middle parameter(s) vary  *)
  548.     filter_[n_, g_List, rest___][ upsample_[l_, args___, n_][x_] ] :>
  549.         upsample[l, args, n][
  550.             resampledfilter[filter, l, n, g, rest][x] ] /;
  551.         isLSIfilterQ[filter] && isUpsamplerQ[upsample] &&
  552.         resampledcoeffs[l, g] ,
  553.  
  554.     downsample_[m_, args___, n_][ filter_[n_, g_List, rest___][x_] ] :>
  555.         resampledfilter[filter, m, n, g, rest][
  556.             downsample[m, args, n][x] ] /;
  557.         isLSIfilterQ[filter] && isDownsamplerQ[downsample] &&
  558.         resampledcoeffs[m, g] ,
  559.  
  560.  
  561.     (*  Polyphase forms of FIR interpolation and decimation structures  *)
  562.     FIR[n_, coeffs_, rest___][ Upsample[l_, n_][x_] ] :>
  563.         PolyphaseUpsample[l, FIR[n, coeffs, rest], n][x] ,
  564.  
  565.     Downsample[m_, n_][ FIR[n_, coeffs_, rest___][x_] ] :>
  566.         PolyphaseDownsample[m, FIR[n, coeffs, rest], n][x]
  567.  
  568. }
  569.  
  570.  
  571. (*  SystemRewriteRules --  62 rules (systems only; no signals allowed)  *)
  572. SystemRewriteRules =
  573.     Join[RulesBasedOnProperties, SingleRateRules, MultirateRules]
  574.  
  575.  
  576.  
  577.  
  578. (*  A P P L Y I N G     T H E     R E W R I T E     R U L E S  *)
  579.  
  580. (*  The Rewrite rules base can be driven by several routines.  *)
  581.  
  582. (*  Recursive Rewrite rules or when all else fails    *)
  583. (*    op_[p__][args__] :> op[p] [ f[args] ] ,        *)
  584. (*    op_[p__] :> op[p] ,                *)
  585.  
  586. Options[ SPRecursiveRewrite ] := { Dialogue -> False }
  587.  
  588. (*  SPRecursiveRewrite   *)
  589. SPRecursiveRewrite[e_, options___] :=
  590.     rewrite[e, Join[ToList[options], Options[SPRecursiveRewrite]]]
  591.  
  592. rewrite[e_, oplist_] := rewrite[e, oplist, Replace[Dialogue, oplist]]
  593.  
  594. rewrite[e_, oplist_, False] := ReplaceRepeated[e, SystemRewriteRules]
  595. rewrite[e_, oplist_, True] := RewriteTracer[e, oplist]
  596. rewrite[e_, oplist_, All] := RewriteTracer[e, oplist]
  597. rewrite[e_, oplist, dialogue_] := Message[Dialogue::notvalid, dialogue]
  598.  
  599.  
  600. (*  RewriteTracer  *)
  601. RewriteTracer[e_, oplist_] :=
  602.     Block [    {difference = True, oldexpr, newexpr},
  603.         newexpr = e;
  604.         Print[e];
  605.         While [ difference,
  606.             oldexpr = newexpr;
  607.             newexpr = Replace[oldexpr, SystemRewriteRules];
  608.             difference = ! SameQ[oldexpr, newexpr];
  609.             If [ difference,
  610.                  Print["which becomes"]; Print[newexpr] ] ];
  611.         newexpr ]
  612.  
  613.  
  614. (*  E N D     P A C K A G E  *)
  615.  
  616. End[]
  617. EndPackage[]
  618.  
  619. If [ TrueQ[ $VersionNumber >= 2.0 ],
  620.      On[ General::spell ];
  621.      On[ General::spell1 ] ];
  622.  
  623.  
  624. (*  W R I T E     P R O T E C T I O N  *)
  625.  
  626. Block [    {newfuns},
  627.     newfuns = { SPRecursiveRewrite };
  628.     Combine[ SPfunctions, newfuns ];
  629.     Apply[ Protect, newfuns ];
  630.     Protect[ SystemRewriteRules ] ]
  631.  
  632.  
  633. (*  E N D I N G     M E S S A G E  *)
  634.  
  635. Print["System rewrite rules have been loaded."]
  636. Null
  637.